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Identifying and understanding interacting systems that can host non-Abelian topological phases with frac¬ 
tionalized quasi-particles have attracted intense attention in the past twenty years. Theoretically, it is possible 
to realize a rich variety of such states by coupling two Abelian fractional quantum Hall (FQH) states together 
through gapping out part of the low energy degrees of freedom. So far, there are some indications, but no robust 
examples have been established in bilayer systems for realizing the non-Abelian state in the past. Here, we 
present a phase diagram of a double-layer bosonic FQH system based on the exact diagonalization and density- 
matrix renormalization group (DMRG) calculations, which demonstrates a potential regime with the emergence 
of the non-Abelian bosonic Moore-Read state. Starting from the Abelian phase with four fold topological de¬ 
generacy at weak coupling, with the increase of interlayer tunneling, we find an intermediate regime with a 
three-fold groundstate degeneracy and a finite drag Hall conductance. We find different topological sectors in 
consistent with Moore-Read state by inserting different fluxes in adiabatic DMRG study. We also extract the 
modular S— matrix, which supports the emergence of the non-Abelian Ising anyon quasiparticle in this system. 

PACS numbers: 73.43.Cd,03.65.Ud,05.30.Pr 


I. INTRODUCTION 

The topological states of matter with emergent fractional¬ 
ized quasiparticles have attracted the intense attentions in the 
past two decades'’^. The statistics of fractionalized quasi¬ 
particles fall into two broad categories: Abelian^ and non- 
Abelian"*^. Interchanging two Abelian quasiparticles changes 
the groundstates by a non-trivial phase factor, whereas inter¬ 
changing two non-Abelian quasiparticles rotates the ground- 
states within a set of degenerated groundstate manifold and 
the final state will depend on the order of operations being 
carried out. Since the quasipartices have unique characteri¬ 
zation to a given topological order, it is possible to identify a 
topological ordered state through the properties of quasiparti¬ 
cles. In particular, compared to the Abelian quasipartcles, the 
non-Abelian quasiparticles are fundamentally important for 
our understanding of the emerging physics, which also have 
potential application for the fault-tolerant topological quan¬ 
tum computation^’®. Currently, the most promising platform 
to investigate the fractional statistics is the fractional quantum 
Hall (FQH) systems, where most of the observed FQH states 
carry Abelian quasiparticles. The prominent candidates which 
may host non-Abelian quasiparticles include the Moore-Read 
(MR) Pfaffian state'^ at the filling factor zz = 5/2^ and Read- 
Rezayi (RR) state® at i> = 12/5'®, which are the FQH states 
for electron systems subject to strong magnetic field. 

FQH states are also observed in the double-layer 
systems"’'^, which can be described in terms of the two- 
component Halperin states'®’'"'. Interestingly, the non-Abelian 
FQH states may also be realized in the double-layer FQH sys¬ 
tems through tuning interlayer tunneling and interactions'® . 
In such double-layer systems, the Halperin wavefunction upon 
symmetrization over the layer index indeed shows the char¬ 
acteristic features predicted for the non-Abelian states, in¬ 
cluding the counting of edge excitations®' and the quasihole 
states'®. In the theoretical considerations, the non-Abelian 
state may be induced by increasing the interlayer coupling, 


which can gap out the low energy degrees of freedom that are 
antisymmetric about the layer inversionHowever, 
there are no strong evidence that this mechanism has been 
realized in physical systems. The numerical studies based 
on exact diagonalization (ED) on the v = 1/2 FQH state 
for electron systems in double-layer find a large wavefunc¬ 
tion overlap between the ground state and the non-Abelian 
state®®^®®, which is consistent with the symmetrization mech¬ 
anism of the Halperin wavefunction. However, the obtained 
energy spectrum has redundant low-energy excitations with¬ 
out a robust energy gap or the groundstate degeneracy on torus 
geometry®®, which are not consistent with a non-Abelian QHE 
state. Very recently, the double-layer bosonic systems on a 
lattice model with topological flat bands have been studied 
variationally based on the parton construction and Gutzwiller 
projected wavefunction. The non-Abelian MR state has been 
identified by using the topological spin and chiral central 
charge®®. However, it remains an open question whether the 
non-Abelian state is indeed the ground state of the micro¬ 
scopic Hamiltonian. As the possible MR state in this double¬ 
layer system appears to be very weak, the accurate simulations 
for large systems and the systematic numerical characteriza¬ 
tion of the topological features are highly desired to pin down 
the nature of the intermediate phase in these systems. 

In this paper, we study a bosonic double-layer system®® 
with each layer in the v = 1/2 Laughlin state in the decoupled 
limit using the ED and density-matrix renormalization group 
(DMRG) calculations. With tuning the tunneling and inter¬ 
action, we find that the four-fold degenerate ground states in 
the decoupled limit evolve to the three lowest-energy states, 
which are symmetry states with respect to the layer inver¬ 
sion. These low energy states are separated from the higher 
energy spectrum by a finite gap in an intermediate parameter 
regime. To identify the nature of the intermediate phase, we 
design and perform different flux insertion simulations®®’®"', 
which can identify total Hall and drag Hall conductances®®. 
We And that the intermediate phase is characterized by the 
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quantized charge Hall conductance and non-zero drag Hall 
conductance, which is distinguished from the Abelian phase 
with the quantized charge Hall conductance and zero drag 
Hall conductance. The flux insertion studies indicate that 
with growing interlayer tunneling, the system evolves from 
the two-component Abelian phase to a one-component phase 
with the nonzero drag Hall conductance, which is consistent 
with a non-Abelian state. Furthermore, we calculate the mod¬ 
ular iS-matrix using ED on finite-size clusters and the obtained 
5-matrix fully agrees with that of the MR state. 

The remaining of the paper is organized as following: In 
Sec. II, we introduce the double-layer lattice model built 
from the topological flat-band (TFB) models. In Sec. Ill, 
we present our phase diagram determined by the energy spec¬ 
trum, charge Hall conductance (Chern number) and drag Hall 
conductance. In Sec. IV, we present the details of numer¬ 
ical results for the quantum phase diagram. We show the 
evolution of the energy spectrum with tunneling obtained 
from ED calculation. The robust groundstate degeneracy is 
checked by inserting charge and spin fluxes in ED. We further 
establish the quantum phase diagram based on the newly de¬ 
veloped large scale adiabatical DMRG simulations, where we 
identify the quantized Chern number and finite drag Hall con¬ 
ductance by inserting the fluxes in the possible non-Abelian 
region. Furthermore, the modular matrix is calculated to sup¬ 
port the emergence of non-Abelian Ising anyon in the possi¬ 
ble non-Abelian regime. In Sec. V, we discuss the topological 
trivial phase in our phase diagram. Finally, in Sec. VI, we 
summarize our main results and discuss open questions. 


II. THEORETICAL MODEL 


We consider a double-layer system composed from two sin¬ 
gle layer topological flat band (TFB) models^®^^, which can 
be generally written as: 

n = n^ + ni + nt + nu, (i) 

= ~ + H.c. , 

(rr') 

"Ht = ^ + H.C. , 

r 

Hu = U± rirlUrf, 

r 

where H^^) denotes the hopping terms in the top layer (bot¬ 
tom layer), Ht describes the interlayer tunneling and Hu 
is the interlayer interaction. s(^r,s) (s =t or |) creates 
(annihilates) a hard-core boson at site r. We consider the 
TFB model on the square lattice and select the phase fac¬ 
tor (pr'r corresponding to half flux quanta per plaquette^^’^^. 
The intralayer hopping terms in include the nearest- 

neighbor coupling J(rr') = 1-0 (energy scale), the next- 
nearest-neighbor coupling = 0.2941, and the next- 

next-nearest-neighbor coupling = —0.2061, which 

give a TFB in each layer with the flatness ratio around 28. 
Such a model can realize the FQH effect for hardcore bosons 
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FIG. 1: (Color online) Quantum phase diagram of the double-layer 
system for the Hamiltonian Eq. (1) with changing the interlayer tun¬ 
neling t± and interaction U±. The blue color filled region represents 
the topological Abelian phase with the fourfold degenerate ground 
states and the topological Chern Number = 1. The green col¬ 
ored squares indicate an intermediate topological order phase with 
threefold degeneracy and C"^ = 1, which is identified as the one- 
component non-Abelian Moore-Read state by the finite drag Hall 
conductance and the modular S-matrix. The yellow dots represent 
a solid phase that breaks lattice symmetries. 


or interacting fermions without a magnetic field due to the 
nontrivial Chern number (C = 1) carried by the topological 
band and the reduced kinetic energy^*^^. 

We consider a finite size system with 2 x Nx x Ny sites 
(Nx X Ny sites for each layer) and the total filling factor of 
the lower TFB v = Np/Ns = {N^ + N^)/Ns = 
where Np is the total number of hardcore bosons and Ng is 
the number of single-particle states of the lower TFB. IN the 
following, we will refer the layer index as the (pseudo)spin in¬ 
dex for convenience. In the absence of tunneling, the system 
reduces to the decoupled two layers with each layer at the 1 /2 
filling. While the ED calculations on torus geometry are lim¬ 
ited to the system with 40 sites, DMRG"^® allows us to study 
the systems up to 256 sites on the cylinder which have the pe¬ 
riodic boundary in one direction and the open boundary in the 
other direction. We keep up to M = 3600 states in DMRG 
calculations, which give accurate results. 


III. PHASE DIAGRAM 


Our main results are summarized as the phase diagram in 
the t± — U± plane as shown in Fig. 1. We find that the Abelian 
FQH state of single layer {t± = U± = 0) is stable against the 
weak interlayer tunneling. In the double-layer system, this 
Abelian phase is characterized by the robust fourfold degen¬ 
erate ground states, the quantized total charge Chern num¬ 
ber = 1, and the vanishing-small drag Hall conductance. 
With the increase of interlayer tunneling on finite-size sys¬ 
tem, the energy of the single antisymmetric ground state splits 
from the other three symmetrized ground states. The higher 
energy spectrum has a gap from the threefold symmetrized 
ground states in an intermediate region (for = 0, 
0.2 < < 0.3). In both the Abelian phase and the interme- 
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FIG. 2: (Color online) Energy spectrum evolution with the interlayer 
tunneling tx on 2 x 4 x 4 torus geometry using ED calculation for 
(a) U±_ = 0.0 and (b) U±_ = 0.5. Among the four near degenerating 
ground states, the three symmetric (S1,S2,S3) states are labeled by 
the black, green, and red circles, while the anti-symmetric (AS) state 
is labeled by the blue stars. 


diate region, we find that the charge Chern Number is always 
quantized as (7° = 1; however, the drag Hall conductance 
jumps from the vanishing-small value in the Abelian phase to 
the almost saturated value in the intermediate region, indicat¬ 
ing that the system evolves from a two-component to a one- 
component QHE state. The phase boundary between Abelian 
and the possible non-Abelian phase is determined by drag Hall 
conductance reaching half of its saturated value (~ 0.125) 
based on DMRG simulation. Interestingly, using the threefold 
ground states, we extract the modular 5—matrix, which sup¬ 
ports a non-Abelian MR state with the emergent Ising anyon 
quasiparticle and the corresponding fusion rule. In the larger 
tunneling regime {t± > 0.3), the system becomes a topologi¬ 
cally trivial state with vanishing Chern numbers (C"^ = 0). In 
the large U_l region, we find a charge density wave state that 
breaks lattice symmetries. 


IV. TOPOLOGICAL NONTRIVIAL QUANTUM HALL 
PHASES 

A. Energy spectrum on torus 

We first use ED to study the evolution of energy spectrum 
with tunneling fx at (7x = 0 on the 2x4x4 torus system 
At t± = 0, both layers have the bosonic v = 1/2 Laughlin 
states with twofold degenerate ground states. Thus, the en¬ 
ergy spectrum of the whole system has a fourfold degeneracy 
separated from higher energy levels by a robust gap. As the 
double layer system has a layer inversion symmetry 
the fourfold ground states can be divided into two groups: the 
symmetric group with three states Esi, Es 2 , Ess shown as 
open circles and the anti-symmetric group with a single state 
Eas represented by blue stars as shown in Fig. 2. By in¬ 
creasing the interlayer tunneling, the groundstate degeneracy 
is lifted gradually. The energy of the anti-symmetric state 
grows rapidly and merges into the high energy continuum 
at t± Ki 0.25. With further increasing fx. two of the three 
symmetric states also merge into high energy continuum at 
t± Ri 0.36. In the intermediate regime 0.20 0.36, the 

three symmetric states have close energies separating from the 


higher energy anti-symmetric state by a finite gap. Although 
there is a finite splitting between the lowest three states due to 
the finite size effect allowing topological distinct states being 
coupled together, the low-energy spectrum of the intermediate 
region implies a possible threefold degenerate ground states 
protected by a gap in the thermodynamic limit. In the pres¬ 
ence of interlayer interaction, we observe the similar results 
as shown in Fig. 2(b) for U± =0.5. Thus, the phase region 
with tx ~ 0.25, which has the maximum finite-size gap, may 
be the suitable phase regime for observing the possible non- 
Abelian Moore-Read state. 


B. Elux insertion based on ED 

We check whether the threefold degeneracy is robust by 
considering the response of the low-energy spectrum to the 
flux insertion. To induce the flux, we impose a twisted 
boundary condition in the y-direction: (ri + Nyy\'^g^) = 
(ri|e®®!/'^o(3) where is the many-body state with 

boundary phase 0y and Pauli matrix a^^s) ^^ts on the layer 
degrees of freedom of the particle at the position r^. The 
twisted boundary is equivalent to threading a flux in the hole 
of a torus along the i-direction"^^. In the double-layer system, 
we introduce two kinds of boundary conditions: charge flux 
(tTo) and spin flux (0-3)^^’'*^’'^®. In the charge and spin flux, the 
boundary phases in the top and bottom layers have the same 
(0^ = 0^) and the opposite (0^ = —0y) signs, respectively. 
For the topological states, the degenerate ground states should 
remain gapped without level crossing with the higher energy 
levels in the charge flux insertion. Therefore, the charge flux 
insertion can be used to identify the near degenerate ground 
states from the low-energy levels. Furthermore, it is expected 
that a two-component double-layer system will have similar 
responses to the charge and spin fluxes, while a coupled one- 
component system has the different responses. 




FIG. 3: (Color online) Energy spectrum evolution with (a-b) insert¬ 
ing a charge flux 0y = 9y and (c-d) a spin flux 6/ = —9/ for 
(fx,Ux) = (0.1,0.5) and (fx, Ux) = (0.25,0.5). The lowest 
five energy levels are shown. 
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(a) t=0.00,U=0.5,charge flux (b) t=0.25,U=0.5,charge flux 
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FIG. 4: (Color online) Real-space configuration of charge accu¬ 
mulation (nj + ri\) after adiabatically inserting a unit charge flux 
9] ^ = 2-k for (a) (f_L, Ui_) = (0.00, 0.5) and (b) (0.25, 0.5). 

The area of the circle is proportional to the accumulation ampli¬ 
tude and the red (blue) color represents the positive (negative) value. 
The results are calculated on 2 x 8 x 16 cylinder using DMRG 
by keeping 3600 states, (c) Real-space configuration of charge ac¬ 
cumulation {n\) after adiabatically inserting a quantized drag flux 
9^ — 2'K{9y = 0). (d) The Charge Hall conductance (red squares) 
and drag Hall conductance (black circles) versus tunneling t±. The 
results for drag flux are calculated on 2 x 8 x 12 cylinder using 
DMRG by keeping 3600 states. 


Fig. 3 shows the ED results of the evolution of low-energy 
spectra with inserting charge and spin fluxes for the system 
with weak and intermediate tunnelings. For = 0.10 and 
C/_L =0.5, the four lowest-energy states are always protected 
by a gap with tuning the charge (Fig. 3(a)) or the spin flux 
(Fig. 3(c)), which are consistent with a two-component topo¬ 
logical state with fourfold degeneracy. In the possible non- 
Abelian phase for (fj^, C/_l) = (0.25,0.5), the lowest three en¬ 
ergy levels are always separated from the higher energy spec¬ 
trum by a small gap in the charge flux insertion (Fig. 3(b)). 
By inserting a spin flux (Fig. 3(d)), we find that the lowest 
three levels connect with higher energy states. These obser¬ 
vations indicate that the intermediate region is a strongly cou¬ 
pled one-component topological state with threefold degener¬ 
acy. However, the nature of the intermediate regime can not be 
established based on these ED calculations. As shown in the 
phase diagram Fig. 1, we determine the phase boundary of the 
Abelian phase where the fourfold degeneracy of energy spec¬ 
trum is being destructed, and a possible non-Abelian phase 
may emerge. We will present the full evidence of the nature 
of the non-Abelian state below. 


C. Flux insertion based on DMRG 

To further investigate the topological properties of the 
double-layer system on larger size, we use DMRG to study 
the cylinder system by adiabatically threading a charge or drag 
flux. The drag flux is realized by introducing the twist bound¬ 
ary phase in just one layer, which can induce a Hall response 


in its own layer and also drag the particles in the other layer. 
Theoretically, the drag Hall conductance and its connection 
to the topological Chern number matrix^^’^^’"^* has been es¬ 
tablished before. Conventionally, one obtains such topologi¬ 
cal Chern invariants based on ED calculations^^’"^^. Very re¬ 
cently, the flux insertion has been introduced in the large- 
scale DMRG simulation on cylinder systems^'^^^, which can 
be used to detect different Hall conductances. 

Very interestingly, the inserting fluxes method we establish 
here for double layer systems can also be used to access all the 
topological sectors in the system. This argument can be easily 
understood in the decoupled limit. Starting from the ground 
state without any flux, we insert the charge flux by adiabati¬ 
cally increasing the twist boundary phase in the closed bound¬ 
ary along the y axis = 0^ = 0 —>■ 27r. By inserting 27r flux, 
the ground state of each layer evolves into a new topologi¬ 
cal sector with a fractional 1/2 charged quasi-particle being 
pumped from one edge of cylinder to the other one^^. Then, 
by adding the drag flux in either layers separately, the system 
would evolve to the other two sectors, which has one more 
pumped charge 1/2 quasi-particle in the layer with drag flux. 
Therefore, we obtain the four topological degenerate ground 
states in the Abelian phase. Qualitatively, this picture in the 
decoupled limit applies to the whole Abelian phase as long as 
the drag Hall conductance is vanishing small, i.e., one layer 
cannot effectively drag the other layer. When the system has a 
transition from two components to one component^*^, the drag 
Hall conductance would jump to a finite value close to the 
saturated value. In the coupled one-component system, the 
drag flux applied to either one of the layer will evolve the sys¬ 
tem to the same topological sector by pumping one fractional 
charged quasiparticle, thus we only obtain three topological 
sectors for coupled phase. 

As shown in Fig. 4(a), after the insertion of one charge flux 
quantum (0j = 0^ = 0 —>■ 27r), a net quantized charge (boson 
number) accumulates at the left edge (nJ + nj. ^ 0.999 is the 
net charge accumulation). The charge accumulation is equiv¬ 
alent to a net charge transfer from the right edge to the left 
edge. According to the fundamental correspondence between 
edge transfer and bulk Chern number^"*, we find a quantized 
charge Chern number = 1 for this Abelian FQH state. 

The quantized charge transfer also persists in the possible 
non-Abelian parameter region. As shown in Fig. 4(b), the net 
charge transfer corresponds to a quantized charge Chern num¬ 
ber C‘^ = 1, which indicates the topological nontrivial nature 
of the possible non-Abelian phase. Nevertheless, the charge 
flux cannot distinguish the Abelian phase from the possible 
non-Abelian phase since both of them share the same = 1. 
To further identify the phase transition between the Abelian 
phase and the possible non-Abelian phase, we consider the 
effects of the drag flux. 

By threading a drag flux quantum in the top layer {9j = 
0 ^ 27r, ej = 0), we observe the particle accumulations 
in the bottom layer^^’^° as shown in Fig. 4(c)). By calcu¬ 
lating the drag Hall conductance as a function of t± as dis¬ 
played in Fig. 4(d), we find a strong enhancement of drag Hall 
conductance at t± ~ 0.15, which coincides with the phase 
boundary of the disappearance of the Abelian phase identified 
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FIG. 5: (Color online) The entanglement spectrum for fx = 
0.25, Ui_ = 0.0 on the Ly — 8 cylinder with (a) no flux, (b) a 
charge flux quanta, and (c) a drag flux (flux imposed on one layer 
only). The entanglement spectrum are grouped by the relative bo¬ 
son number AN of left half cylinder. The results are obtained using 
DMRG by keeping up to 1600 states. 


from the ED energy spectrum calculations. Within the regime 
0.25 < t± < 0.30, the drag Hall conductance approaches 
the saturated value 0.25, which is consistent with an effec¬ 
tive one-component system. Based on the above results, we 
determine the Abelian phase (blue color), the one-component 
non-Abelian phase (green squares) as shown in Fig. 1 . The 
topological trivial solid phase with zero charge Chern number 
is also shown. 

To further identify the possible non-Abelian state, we in¬ 
vestigate the entanglement spectra of each topological sector. 
Physically, the different topological ground states on a cylin¬ 
der are expected to have the different well-defined anyonic 
flux through the cylinder. Thus, the cylinder system with the 
charge flux or the drag flux corresponds to the other two topo¬ 
logically distinct ground states besides the vacuum state. To 
explicitly demonstrate the ground states with different any¬ 
onic flux on cylinder geometry, we bipartite the cylinder into 
two halves, and observe entanglement spectrum^'^ to distin¬ 
guish the different topological sectors. As shown in Fig. 5, 
we show the entanglement spectrum for the vacuum ground 
state in Fig. 5(a), the new ground state obtained by insert¬ 
ing a charge flux quantum in Fig. 5(b), and the ground state 
obtained by inserting a drag flux in Fig. 5(c). These three 
ground states are anticipated to have one-to-one correspon¬ 
dence with identity, fermion, and Ising anyon sectors, respec¬ 
tively. We also calculate the momentum dependence of the en¬ 
tanglement spectra in each (7(1) quantum number sector with 
different relative boson number AN, and obtain the counting 
of the leading eigenvalues in the entanglement spectra®^. The 
obtained results are similar to those of coupled two Faughlin 
1 / = 1/2 states. Due to the calculation limit, = 8 (16 lat¬ 
tice sites in the x direction) is the largest width we can reach 
convergence in our DMRG calculations, which gives four mo¬ 
mentum quantum numbers K = 0, 7 r/ 2 , tt, 37 r /2 in each AN 
sector. Although we observe that a very small entanglement 
gap opens up in AT = tt and 37 r /2 sectors between the ex¬ 
pected counting for non-Abelian MR state and the other part 
of entanglement spectrum, we cannot determine if the gap will 
survive in the thermodynamic limit or it is a finite size effect. 


Since all other results support the non-Abelian QHE state, we 
believe this result is due to the finite size effect and we leave 
this part for the future study. 


D. Modular matrix 

From the above observation of DMRG, we And that the in¬ 
termediate phase region appears as a one-component topologi¬ 
cal nontrivial phase. Here we calculate the modular 5—matrix 
using the near degenerate threefold states in the ED energy 
spectrum to further investigate the nature of the possible 
non-Abelian phase*. Modular 5—matrix encodes the infor¬ 
mation of the quasiparticle statistics including quantum di¬ 
mension and fusion rules^®“^^, which has been successfully 
used to identify various Abelian and non-Abelian topological 
orders^®^®^. To calculate the modular matrix, we follow the 
method based on the minimal entangled states (MESs)®°. The 
MESs are the eigenstates of the Wilson loop operators with a 
definite type of quasiparticle^®. Thus, the modular transfor¬ 
mations on the MESs give rise to the modular matrix®**. 

Here we show the results at t± = 0.25, U± = 0.0 as an ex¬ 
ample. We denote the three lowest-energy states in ED spec¬ 
trum as (/ = 1 , 2 ,3), from which we can form the general 
superposition states as, 

l'I*(c.,C3,0..03)) = Cl 16) + C 26 ^= 16) + C36^^ 16) 

where ci, C 2 , C 3 , (j) 2 , />3 are real superposition parameters. For 
each state Ifll), we construct the reduced density matrix and 
obtain the corresponding entanglement entropy. To And the 
MESs, we optimize the superposition parameters to And the 
minimum entanglement entropy. As shown in Fig. 6 (a), we 
show the entropy profile of Ifll) with the optimized parameters 
{(j) 2 , 'PV) middle cut along the i-direction. We And the 

first global MES |S{) with the entropy S ^ 3.10 at the posi¬ 
tion pointed by the red arrow. The second MES (blue ar¬ 
row) and the third MES jSg) (green arrow) can be determined 
in the state space orthogonal to |S{), as shown in Fig. 6 (b). 
Finally, we can obtain the modular matrix S = ex- 
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FIG. 6 : (Color online) Left: Contour plot of entanglement entropy of 
I^'ci,c 2 , 23 ,^ 2 . 03 ) on 2 X 4 X 4 system. We show entropy profile ver¬ 
sus C2, C3 (ci = a/I ~ Cj — C 3 ) by setting optimized 4>2,'p3- Three 
nearly orthogonal MESs are marked by red, green and blue arrows 
in surface plot. The cyan dashed line represents the states orthogonal 
to the first MES (red arrow). Right: Entropy for the states along the 
cyan dashed line as shown in left figure. 
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FIG. 7: (Color online) (a) Entanglement spectrum for t± = 
0.08, U± = 1.2. (b) The charge density distribution of (nj) (top 
layer) and (nj) (bottom layer) in real-space. Here we only show two 
columns in the bulk in each layer. The two-column pattern is periodic 
along the cylinder direction. The area of the circle is proportional to 
the amplitude of charge density. The results are calculated on the 
2 X 8 X 16 cylinder using DMRG with keeping 1600 states. 


V. TOPOLOGICAL TRIVIAL PHASES 


Besides the two topological non-trivial phases in the phase 
diagram Fig. 1, we also find the topological trivial phases 
at the large U± parameter region. Here we show that there 
is a charge density wave phase at large U± region, as shown 
by the yellow triangular in Fig. 1. As shown in Fig. 7 (a), 
the entanglement spectrum has a large weight on the lowest 
eigenvalue, which is a feature of a solid phase. In Fig. 7(b), 
we demonstrate the charge density wave pattern in real space. 
Along the cylinder axis (a;—direction), the unit cell (enclos¬ 
ing 4 sites) is doubled in each layer due to the charge density 
pattern. In the y—direction, the charge density pattern has a 
period of two. At each site, the top and bottom layers have the 
opposite density pattern, which leaves (nj + n^) a constant. 


tracted from the overlap between the MESs for two noncon- 
tractible partition cut directions^®: 
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where represents the theoretical prediction from the 
SU{2)2 Chern-Simons theory^^^^®. determines the 

quasiparticle quantum dimension as = 1, = 1, 

da = and non-trivial fusion rule as tr x cr = i + ip, 
where 1 represents the identity particle, ip the fermion-type 
quasiparticle, a the Ising anyon quasiparticle. Thus, the nu¬ 
merical extracted modular 5—matrix identifies the interme¬ 
diate topological phase with threefold ground state degener¬ 
acy as the non-Abelian MR state with the emergence of the 
Ising anyon quasiparticles satisfying the non-Abelian fusion 
rule (533 ~ 0)«-«. 

Generally speaking, to uniquely determine a topological or¬ 
der, one needs both the modular S and U matrices®®. From 
modular U matrix, one can access the chiral central charge 
c and topological spin of each quasiparticle, which distin¬ 
guish the non-Abelian MR state from double-layer Abelian 
Laughlin states. For example, the chiral central charge of non- 
Abelian MR state is c = 3/2, while the double-layer Laughlin 
state has c = 2. Unfortunately, in the current lattice model 
(in Eq. 1 ), the MES route can not give the S and U matrix 
together since there is no rotation tt/S symmetry here®®. Re¬ 
cently, we note that it has been proposed a general method, 
named momentum polarization®^, to extract the quasiparticle 
statistics in modular U matrix. In this method, one needs to 
perform a finite-size scaling on Ly, then the statistics informa¬ 
tion in modular Z/f matrix can be extracted from Ly ^ 0 limit. 
Unfortunately, in current DMRG calculation we are limited to 
the system sizes Ly = 4, 6,8 due to the computational capa¬ 
bility. Thus a reliable application of momentum polarization 
method here is very challenging, which we leave for the future 
study. 


VI. SUMMARY AND DISCUSSION 


We have studied a bosonic double-layer system on a square 
lattice using ED and DMRG calculations. Through the 
studies of the energy spectrum, the flux insertion on cylin¬ 
der, and the modular matrix, we find numerical evidences 
for a non-Abelian Moore-Read state emerging from the bi¬ 
layer Halperin states through gapping out the interlayer anti¬ 
symmetric state. Although this practically powerful route to 
a variety of non-Abelian quantum states has been introduced 
theoretically for decades based on parton construction and 
field theory^®^®®’®®^®®, there were limited numerical evidence 
to support the realization of non-Abelian state in microscopic 
systems. Our numerical calculations rely on the insertion of 
charge and drag fluxes, which allow us to detect the quan¬ 
tum phase transition from a two-component topological state 
to a one-component state characterized by the onset of the fi¬ 
nite drag Hall conductance. In combining with the modular 
matrix simulation for the quasiparticle statistics, we identify 
the nature of the intermediate t± (with threefold near degen¬ 
eracy) phase as the non-Abelian Moore-Read state, although 
this state is relatively weak and the entanglement spectrum 
does not show a robust entanglement spectrum gap for the 
counting associated with the non-Abelian state. 

We have also explored the nature of the quantum phase 
transition between the Abelian phase and the possible 
non-Abelian phase. We have studied quantities such as entan¬ 
glement entropy and the wavefunction overlap. We find all 
of these quantities of the groundstate change smoothly when 
the system crosses the phase boundary. This indicates that 
the phase transition is either weakly first order or continuous 
transition. Moreover, it would be particularly interesting to 
study the possibility of realizing the non-Abelian phase from 
coupled bilayer Halperin states in fermionic systems, which 
will be investigated in the future work. 

Note added. Upon finalizing the manuscript we noticed 
several preprints focusing on double-layer v = 1/3 + 1/3 
fermionic systems®®^®®. 
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